****************************************************************************
**		Authors: 	Or Tuttnauer & Liran Harsgor
**		Purpose: 	Produce Figure 6 in Kedar, Harsgor & Tuttnauer (JOP)  
**		input:		KHT_countrylevel.dta
*****************************************************************************

log using "figure 6.log", replace

set more off
*Import presidential candidates data
use "presidential_data.dta", clear

*Generate two vars for C+G analysis
gen share_sq=candidate_share^2
bysort cyear: egen sum_squares=total(share_sq)
gen enpres=1/sum_squares
gen proximity = 2 * abs(((legis_year-prev_pres)/(next_pres-prev_pres))-.5) 

*reduce to one obs per country and merge the two generated vars into our dataset
egen country_u = tag(cyear)
keep if country_u==1
keep cyear enpres proximity
save KHT_temp.dta, replace

use "KHT_countrylevel.dta", clear
merge 1:1 cyear using KHT_temp.dta
drop if _merge==2
drop _merge

recode enpres proximity (.=0)

** Analysis **


* Interaction terms *
gen elavgdm = elffrn*lnavgdm
gen eupper   = elffrn*upper
gen esr      = elffrn*sr
gen presprox = enpres*prox


*** 10% and 90% for center and dispersion variables
sum sr, detail
scalar sr10=r(p10)
scalar sr90=r(p90)

** Generate for scatter
	gen disp=0.3
	gen dispavgdm=lnavgdm+(runiform()-.5)/10
	gen pipe="|"

	

reg enpv elffrn lnavgdm sr elavgdm esr upper enpres prox eupper presprox

* Save coefficients and V-CV matrix *	
matrix b=e(b)
matrix V=e(V)
scalar b1=b[1,1]  		/* first coefficient (elffrn) */
scalar b4=b[1,4]		/* first interaction coefficient */
scalar b5=b[1,5]		/* second interaction coefficient */
scalar varb1=V[1,1]		/* variance of b1 */
scalar varb4=V[4,4]		/* variance of b4 */
scalar varb5=V[5,5]		/* variance of b5 */
scalar covb1b4=V[1,4]	/* covariance of b1 and b4 */
scalar covb1b5=V[1,5]	/* covariance of b1 and b5 */
scalar covb4b5=V[4,5]	/* covariance of b4 and b5 */

generate MVZ= ((_n-1)/14)

gen conbx= b1+b3*MVZ
gen consx= sqrt(varb1+varb3*(MVZ^2)+2*covb1b3*MVZ)
gen a= 1.96*consx
gen upperx= conbx+a
gen lowerx= conbx-a

	
gen conbx10= b1+b4*MVZ+b5*sr10
gen consx10= sqrt(varb1+varb4*(MVZ^2)+2*covb1b4*MVZ+(sr10^2)*varb5+ ///
2*sr10*covb1b5+2*MVZ*sr10*covb4b5)
gen conbx90= b1+b4*MVZ+b5*sr90
gen consx90= sqrt(varb1+varb4*(MVZ^2)+2*covb1b4*MVZ+(sr90^2)*varb5+ ///
2*sr90*covb1b5+2*MVZ*sr90*covb4b5)

* Confidence intervals
gen a10= 1.96*consx10
gen upperx10= conbx10+a10
gen lowerx10= conbx10-a10
	
gen a90= 1.96*consx90
gen upperx90= conbx90+a90
gen lowerx90= conbx90-a90

* Generate graph
	graph twoway scatter disp dispavgdm if smd==0, ///
						ms(none) mlabel(pipe) mlabposition(0) mcolor(black) yline(0) ///
		||	 line conbx10 MVZ, clpattern(longdash) clwidth(medium) clcolor(black) ///
        ||   line conbx90 MVZ, clpattern(solid) clwidth(medium) clcolor(black) ///
        ||   line upperx10 MVZ, clpattern(longdash) clwidth(thin) clcolor(gray) ///
        ||   line lowerx10 MVZ, clpattern(longdash) clwidth(thin) clcolor(gray)  ///
		||   line upperx90 MVZ, clpattern(sold) clwidth(thin) clcolor(gray)  ///
        ||   line lowerx90 MVZ, clpattern(solid) clwidth(thin) clcolor(gray) ///
			, legend(off) graphregion(fcolor(white) lcolor(white))  ///
            xtitle("Average district magnitude (logged)", size(3)) ///
			ytitle("Effect of ELF on ENPV", size(3))  ///
            scheme(s2mono) xscale(range(0 3.7))	xlabel(0(1)3) ylabel(-4 (2) 14)

graph save "figure 6.gph", replace

log close
